Effect of Shear Flow on the Stabihty of Domains in Two Dimensional 

Phase-Separating Binary Fluids 
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We perform a linear stability analysis of extended domains in phase-separating fluids of equal 
viscosity, in two dimensions. Using the coupled Cahn-Hilliard and Stokes equations, we derive 
analytically the stability eigenvalues for long wavelength fluctuations. In the quiescent state we find 
an unstable varicose mode which corresponds to an instability towards coarsening. This mode is 
stabilized when an external shear flow is imposed on the fluid. The effect of the shear is seen to be 
qualitatively similar to that found in experiments. 
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I. INTRODUCTION 



Phase separating binary fluids form complex patterns 
of domains at late times after a temperature quench into 
an unstable state. The morphology of the domains is de- 
termined by factors such as the volume fractions of the 
two phases, the viscosities of the two phases, and any 
externally applied forces Of particular interest to 

us is the effect of applying an external shear flow to a 
phase separating binary fluid. This question is of tech- 
nological importance because many industrial processes 
involve binary mixtures in a flow field. The final material 
properties depend on the domain morphology, which can 
be strongly affected by the fluid flow. 

At late times after a temperature quench into the two 
phase region of the phase diagram, a phase separating 
fluid consists of domains of the two phases of typical size 
R(t), which coarsen with time generally as a power law 
R{t) oc t" [Qjsj. The presence of a shear flow dramatically 
alters the kinetics of the phase separation. The shear 
flow deforms the domains, interfering with their growth 
so that it competes with the thermodynamic force driving 
the phase separation. Many theoretical 0-0] and exper- 
imental [|-10| studies have investigated the effect of the 
shear flow on the growth of the domains and the exponent 
a. In this work we focus on a different aspect of the effect 
of shear: eventually the binary fluid tends towards a dy- 
namic, nonequilibrium steady state in which the coarsen- 
ing instability is stopped by the shear flow ; 
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The 

morphology in this stationary state is very anisotropic . 
In relatively weak shear, the domains are somewhat de- 
formed, whereas at higher shear they can become highly 
elongated along the flow direction. A "string phase" con- 
sisting of macroscopically long cylindrical domains forms 
when the two phases are both percolated [ p"3p4| ]. This 
is surprising, since a long cylinder of fluid at rest would 
normally break up via the Rayleigh instability |15 1^, a 
hydrodynamic instability. The string phase appears to 
be a fairly robust phenomenon, appearing in both criti- 
cal and off-critical polymer mixtures |0| and in critical 



micellar solutions |T^. Thus, the shear flow both opposes 
the thermodynamic instability driving phase separation 
and stabilizes these highly anisotropic domains against 
hydrodynamic instabilities. 

Our goal is to understand these stabilizing effects of 
shear flow. As a first step towards elucidating these ef- 
fects, we consider a strictly two dimensional system. We 
expect the operative physical mechanisms in the two di- 
mensional fluid to be somewhat different than those in 
the three dimensional case, but the mathematical tech- 
niques and physical insights developed here will be of 
use in the future for three dimensional calculations. We 
consider late times after an initial temperature quench 
into the unstable region of the phase diagram, when 
the system is composed of domains of the two phases 
close to their equilibrium concentrations and separated 
by well-defined interfaces. We will, however, retain the 
dynamics of the concentration field in our analysis, so 
that the interfaces between domains have a finite width 
^. We model the fluid using the coupled Cahn-Hilliard 
and Navier-Stokes (for creeping flow) equations as de- 
scribed in Section O. T his is in contrast to the work of 
San Miguel et. al. |]l8||, who did an analysis of the sta- 
bility of domains in two dimensional binary fluids, using 
only the Navier-Stokes equation and treating the inter- 
faces as mathematically sharp. 



In Section [II we linearize our equations for the general 
case of a system with any number of flat interfaces, and 
develop some useful mathematical machinery. In Section 



IV 



we apply our methods to the case of a single interface, 
and reproduce some well-known results. In Section ^ we 
turn to our main focus, the stability of a single domain 
in the form of a strip (in three dimensions, a flat sheet) 
of one phase, immersed in an infinite region of the other 
phase as illustrated in Fig. |l|. We impose a shear flow 
along the x-direction by applying a constant shear stress 
IIq. In this paper we take the viscosity of the two phases 
to be equal, so that the flow field of the unperturbed 
system is linear. There are two linearly independent per- 
turbations of the lamellar domain along the x-axis. In 
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the "zig-zag" mode the two interfaces fluctuate in phase, 
whereas in the "varicose" or "peristaltic" mode they fluc- 
tuate out of phase. We find that in the absence of the 
shear flow the zig-zag mode is stable, whereas the vari- 
cose mode is unstable to long wavelength perturbations. 



We use a tight-binding approximation to include the ef- 
fect of the shear flow. In Section Vl we observe that the 



shear flow mixes the two modes so that above a critical 
shear rate 7c the lamellar domai n is stable. We conclude 
with some discussion in Section VII. 
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FIG. 1. Geometry of a single lamellar domain of phase a 



II. MODEL EQUATIONS 

We consider a simple binary fluid with one scalar or- 
der parameter the difference in concentration between 
the two components. We use the usual Ginzburg-Landau 
form for the coarse-grained free energy of a symmetrical 
mixture 



F[<i>] 



dr 



2 ° 



4^ 



(2.1) 



where Tq and g are positive so that we are below the co- 
existence curve in the two-phase region. Minimizing the 
homogeneous part of F leads to the values of the concen- 
tration in the two bulk phases at equilibrium: 



$ = ± 



±0 



The equation of motion for $ is the Cahn-Hilliard equa- 
tion with a convective coupling of $ to the velocity field 



dt (5$ 



(2.2) 



Here Af is a concentration-independent mobility. Since 
we are interested in the late stages of phase separation, 
we neglect all thermal fluctuations. The equation for the 
velocity is the Navier-Stokes equation for an incompress- 
ible fluid, generalized to include the coupling of the order 
parameter to the velocity field 119] : 



p— + p(u • V)u = 77V^u 



(5F 

V$— - VP. 
i5<i> 



(2.3) 



In this paper the viscosity r\ will be taken to be indepen- 
dent of $; hence there is a single viscosity for the fluid 
independent of the concentration pattern. The pressure 
P is determined by the incompressibility condition 



V • u = 0. 



(2.4) 



We will consider only low Reynold's number flow so the 
convective term (u • V)u in the Navier-Stokes equation 
can be ignored. We will also assume that the fluid is 
sufficiently viscous that the velocity responds instanta- 
neously to slow changes in we can then neglect the 
inertial term dn/dt and the resulting equations describe 
"creeping" or Stokes flow. The term coupling the concen- 
tration to the velocity field in (2_^) leads to a capillary 
force at i nterf aces, where gradients in $ induce fiuid flow. 
Eq. (U - O) are the same as those of "Model H" (with- 



out the thermal noise terms) used to study critical binary 
fluids These equations have been used extensively 

to study phase separation in binary fluids 

The first step in a stability analysis is to derive the 
steady-state solutions to the equations of motion. With 
the geometry of Fig. |l| in mind, we assume that $ and 
u are functions of y only and look for time-independent 
solutions. The Cahn-Hilliard equation ( p. 2D has steady 
state solutions satisfying 



5F_ 



-fCV^<I> - ro<I> -I- .g$^ = M = constant , (2.5) 



where is the exchange chemical potential. Near a single 
interface, we can take /i = and the concentration has 
the usual "kink" solution 



tanh - 



;tanhy/^, 



(2.6) 
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where the width of the interface between the two coexist- 
ing phases is the thermal correlation length ^ = \/2 K /rp. 
For a system of many lamellar domains, Eq. (Ej) gives 
the concentration profile at each interface when the do- 
main size is much larger than ^. We note that there is a 
surface tension associated with the presence of an inter- 
face, which is just the excess free energy per unit area at 
the interface |2lf| : 



a = K 



dy J 3g 3g 



(2.7) 



In the stationary state in shear flow there is no velocity 
in the y-direction. We impose a constant shear stress Uq 
so the stationary velocity satisfies 



Us = jyx 



(2.8) 



where 



7 



Ho 



is the shear rate. 

It is convenient to rewrite our equations in dimension- 
less form by scaling lengths by the correlation length, the 
concentration by its equilibrium magnitude in the bulk 
phases, and time by the natural diffusion time involving 
the mobility M in the Cahn-Hilliard equation. The ve- 
locity is scaled by the correlation length over the diffusion 
time: 




Mii_2Mr^ 



$1' 



K 



P = P 



Mrl^ 2Mro ' 
2K(j)l ■ 



Note that the new dimensionless correlation length is 
^ = 1. In dimensionless form the equations of motion 
are now 



5$ -- 1-, 



1. 



-V^$ - $ + 



= V^u + f -iv^* - $ + 
= V • u. 



(2.9) 



-VP, (2.10) 



(2.11) 



We see that the system is characterized by the dimen- 
sionless parameter, fj: 



V ■ 



Mgr] AMroT] 



K 



2,ai 



(2.12) 



In dimensionless form, the concentration and velocity 
profiles derived above for a single interface parallel to 
the ffow are 



^s{y) = tanhy , 



Us(y) = 72/x. 



(2.13) 



(2.14) 



The dimensionless shear rate 7 = 7idiff is simply the 
product of the shear rate and the diffusion time tdiff — 
^ IMto, and thus represents a second dimensionless pa- 
rameter that characterizes the strength of the shear flow. 



III. STABILITY ANALYSIS 

In this section we will develop an overall strategy to 
examine the stability of any number of lamellar domains. 
We perform a linear stability analysis about the sta- 
tionary states derived above. We begin by considering 
small perturbations about the stationary solutions (we 
will drop the bars over the dimensionless variables in the 
rest of the discussion, except where noted): 



= $ - 



V = u u. 



(3.1) 
(3.2) 



To linear order in the perturbations and v the equa- 
tions of motion become 



1. 



dy 



1 



' dx 



1 



1. 



Ws{y) 



(3.3) 



W+-^(--V^ + M/,(2/))0y-VP, (3.4) 



= V • V. 



(3.5) 



Here Ws is a function of the stationary concentration 
profile: 



Ws{y) 



-l + 3</>2(2/). 



(3.6) 



For a single interface at y = 0, Ws{y) — 2 ~ 3sech y 
so that the nonconstant part of Ws is isolated near the 
interface. 

In this work we are interested in perturbations along 
the flow direction that are perpendicular to the planar 
interfaces. Any such perturbation can be written as a 
sum over Fourier components along the x-direction, so 
we take our perturbations to have the plane-wave forms 



4> = Hyy"' 



ikx — Lot 



V = v(y)e 



ikx—uit 



(3.7) 
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We will consider long wavelength fluctuations for which 
A;^ <C 1. Note that in the following, we take k to be pos- 
itive, so that k represents the magnitude of the wavevec- 
tor. First we consider the hydrodynamic equations. If 
we substitute the expression fo r v given in (|3.7D into the 
equations of motion for v, Eq. (3^) and Eq. ( |3.5|) , we find 
we can solve them exactly in terms of a Green's function. 
First we introduce the stream function 4*, defined by 



dy ' 



dx 



(3.8) 



The incompressibility condition (3_^) is then automati- 
cally satisfied by The two components of the Navier- 
Stokes equation ( |3.4|) can be used to eliminate the pres- 
sure P, leaving a fourth-order ordinary differential equa- 
tion for ^' = i^iy) exp{ikx — tot): 



1 



(3.9) 



Here primes indicate differentiation with respect to y. 
The boundary conditions are that i/; and its derivative 
vanish at infinity. This equation can be formally solved 
using a Green's function, to obtain the ^/-comp one nt of 
V that is needed in the concentration equation (3.3): 



X Qfc2</'(2/') - \i>"{y') + T^.(2/')0(2/')) ■ (3.10) 

This gives Vy in terms of an integral over 

Next substituting (3.7) into the concentration equation 



(3.3) results in an eigenvalue equation for io(k\. 



;(fc)0 = 



" dy 



ikjycf) 



(3.11) 



dy'^ 



where we have used Ms = ^y. A real, positive value of 
Lj{k) indicates stability (damping) of the perturbation. 
Note that this is essentially an integro-differential equa- 
tion in which Vy acts as an integral operator on <f>. 



We cannot solve Eq. (3.11) exactly so an approximate 
method is needed. To deve lop our calculational approach 
we first consider Eq. ( 3.11 ) without the flow terms: 



Ljd)^ — T - A 

2 \dy^ 



2dy^ 



(3.12) 

This equation is applicable to the perturbations of do- 
mains in a binary solid and was used by Langer to 



descr ibe c oarsening mechanisms in binary alloys. Note 
that ( 3.12| ) has the form 



uj(t) = rF4> , 

where we have defined the operators 
1 / d^ 

r : 



2 ydy^ 



F : 



2dy2 



Wsiy). 



(3.13) 



(3.14a) 



(3.14b) 



If (j>n is the set of eigenfunctions of ( 3.13| ) and we define 
a set of "conjugate" functions by 



r(f)„ 



(3.15) 



then one can show that F and F are Hermitian opera- 
tors (although their product is not), as long as the (/)„ 
and (pn obey periodic boundary conditions or vanish at 
infinity. We note that the eigenvalues ujn are real and the 
eigenfunctions and their conjugates are orthogonal: 



4>*ni{y)4>n{y)dy = O for n 7^ m . 



Then for any pair of trial functions (pa and (po obeying the 
same boundary conditions, wc can find an upper bound 
on the lowest eigenvalue w from a variational relation 



< 



(0o,0o) 



(3.16) 



Here the parentheses again indicate inner products. 

To apply Eq. (3.16) we need a good tri al fun ction 00- It 
is easy to determine an exact solution of ( 3.12| ) in the par- 
ticular case when we have a single flat interface present 
and when (/) is a function of y only [k = 0). We note 
that the system is translationally invariant, so that any 
solution that corresponds to a translation of the interface 
by some amount dy is also a solution. Thus \iy ^ y + dy 
we can write 



(ps{y + dy) = 0s + -r-dy + 
dy 



so it must be that 



— — — sech y 
dy 



(3.17) 



is also a solution. It is easy to verify that this is the 
case, with corresponding eigen value cj = 0. This is the 
lowest lying eigenvalue of Eq. ( 3.12 ) for a system with 
a single planar interface and fc = ||2^ . We can use 
the variational principle ( 3.16| ) to calculate the stability 
eigenvalues near this w = translational mode for more 
general situations by assuming a trial function formed by 
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appropriate linear combi natio ns of the single interface so- 
lution To use Eq. ( 3.16 ) we also need to determine 
the conjugate function ^o- By definition the conjugate 
function satisfies 



1 

2 \dy^ 



- e My) = My) 



(3.18) 



We can easily solve for by using a Green's function, 
with boundary conditions that (j)Q and 
ity. We find 



Q vanish at infin- 



dy'p-k\y 



y 



My') 



(3.19) 



The conjugate function is thus obtain ed by substituting 
the desired trial function into Eq. ( [3.19| ). 

To summarize the results of this section, we have lin- 
earized the equations of motion, expressed them para- 
metrically in terms of the wavenumber k, and solved the 
hydrodynamic equati ons f or Vy as an integral over cj). The 



eigenvalue equation ( 3.11 ) can be solved approxi mate ly 
in the absence of the t wo flow terms (i.e. Eq. ( 3.12] )) 
by evaluating Eq. (3.16) using an appropriate trial func- 
tion. The methods used to include the flow terms will be 
explained in the following sections. 



IV. DISPERSION RELATION FOR A SINGLE 
INTERFACE 

As an example of the variational technique, consider 
the dispersion relation of a single flat interface separat- 
ing semi-infinite domains of the two phases. We initially 
negle ct hydrodynamic effects and focus on solving Eq. 
( 3.12 ) for uj{k). For a single interface located at y — 



our trial solution is exactly (f)o = 4>'g = sech y. The re is 
only one term in F(j)o, since (j)o is a solution to (3.12) for 
k = 0: 



1 £ 1 

J^-Z-o = ( - + 2^' + 2 - 3sech"2/ ) (sech^y) 
1 



= ;^^^sech^u . 
2 ^ 



Using (3.19) for the conjugate function (j)Q, we find 



dy'j exp{-k\y - y'|)sech^y' 
k 



IK 



\y-y\-\ ) sech^y' 



= 2 In cosh y -|- 0(k) , 

k 

where we have expanded the exponential for small k (long 
wavelengths). This expansion is not uniform in y, but is 
justified since the integrand is only nonzero for small y' 
and because we will only need (po for small values of y 



in the subsequent analysis. The normalization integral is 
simply 



i(t>o, M 



dy [ ^ — 2 In cosh J/ + 0{k) ) sech^y 



--2(2-21n2)-|-0(fc). 
k 



Next we apply the variational theorem ( |3.16 ) to obtain 



4/fc-2(2-21n2) + 0(A:) 



(4.1) 



^ k' + 0{k*) 
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where we have retained only the lowest order term in k. 
If we rewrite this relation in dimensional units, we find 



uj^\Dk^^ + 0{k^) 



(4.2) 



where D = Afro is a diffusion constant. This result has 
been obtained previously by Jasnow and Zia ]2^ and by 
Shinozaki and Oono |25|. It also agrees to lowest order 
in k with the perturbative calculation by Bettinson and 
Rowlands The eigenvalue is positive so the single 

interface is stable, at least to long wavelength perturba- 
tions. 

The physics here is straightforward. We know that 
outside a curved interface there is a slight excess concen- 
tration, which is given by the Gibbs-Thomson relation 

M 



RAcj) 



(4.3) 



where a is the surface tension, x is the susceptibility, R is 
the radius of curvature of the interface, and Acj) = 2(j>e is 
the miscibility gap. In our case the curvature of the inter- 
face is l/R — Ak^ where A is the amplitude of the small 
perturbation. The susceptibihty x is X^^ = d^i/d<l) = Tq 
in the bulk phase. The excess concentration due to the 
curvature is therefore 



where we have used (2/7) to eliminate a. This excess 
concentration will occur outside regions of positive cur- 
vature, and there will be a corresponding lack of concen- 
tration in regions of negative curvature, creating a con- 
centration gradient along the x-axis. The flux across the 
interface caused by this gradient is roughly vA.(j) where v 
is the velocity of the interface. That velocity, in turn, is 
just the rate of change of the amplitude A of the pertur- 
bation, so 

dA 
dt 
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The concentration gradient is Wcj) kScj); putting every- 
thing together, we find 

so that uj ^ Dk^^ as advertised. 

We can include the lowest order hydrodynamic effects 
on the dispersion relation by performing a perturbative 
calculation to fi rst o rder in k. We write the full eigen- 
value equation (3.11) in the form 



rF(t) + V(j) = uj(l), (4.4) 
where the "unperturbed" problem is simply Eq. ( 3.12| ): 



2 [dy^ 

LUo4>0 



Wsiy) 00 



with luq = fc^/6 and (f>o = sech^y from the variational 
result (note these solutions are exact for k = 0). The 
perturbation V contains the flow terms 



V 



dy 



ikjycj} . 



We expect Vy to be proportional to a power of k (since 
for fc = there should be no induced velocity in the y- 
direction), so V itself is proportional to a power of k and 
is therefore small for long wavelengths. Expanding lu and 
(j) in powers of fc, and multiplying Eq. (4.4) on the left 
by the corresponding left eigenvector (jj, one can show in 
the usual way that the lowest order correction to uj in 
perturbation theory is 



1^(1)0, V(j)o 



(4.5) 



We s olve for the velocity field by substituting (po into 



r dy' {l + k\y-y'\)c~^\y~y'\^ec\v^y' 

4?yfc 



X I -fc^sech^y' 
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dy' (sech-^y' + 0(fc2)) 



(4.6) 



where we have again expanded the exponential for small 
fc. We find that to lowest order Vy is linear in fc, so that 
overall V (x k. Since our first order perturbative result 
will only be good to 0{k), we only need the exact part of 
the result to the unperturbed problem (recall the varia- 
tional result is 0{k^)), for which luq — 0. In the reference 



frame in which Us(y) = j y, the integral over the convec- 



tive term ikjycj) in (4^) vanishes, so that we obtain a 
single term in the first order correction to uj from the Vy 
term: 



UJ 



LUq ■ 



»0, (Po 



= 



0o,(fc/677 + O(fc2)) 00 



brj 



(4.7) 



since (f>'g = sech^y = (/jq for a single interface. If we restore 
the units in this result we obtain 



.= ^-fO(fc^) 



(4.8) 



where a is the surface tension. This is a well-known re- 
sult for the damping of long wavelength capillary waves 
on a planar interface between two liquids, in the limit 
that the viscosity is sufficiently large that inertial effects 
can be neglected ||2^ . 



CALCULATIONAL METHOD FOR A 
LAMELLAR DOMAIN 



We now turn to the stability of a lamellar domain of 
one phase immersed in the other phase, so that we have 
two interfaces in the system as in Fig. H When the spac- 
ing A between the two interfaces is at least a few correla- 
tion lengths (note that we continue to work with scaled 
variables), A 3> 1, the stationary concentration profile is 



tanh(y + A/2) -oo <y <0 
- tanh(j/ - A/2) <y < +00 



(5.1) 



where we have arbitrarily taken the a phase with equi- 
librium concentration (jja — +1 to be in the middle, with 
layer thickness A. In this expression we have set the ex- 
change chemical potential /j, to zero. More accurately, we 
can calculate uas follows. The stationary solution that 
satisfies Eq. ( ^.5[ ) is 

<j)s = tanh(?/ + A/2) - tanh(?; - A/2) + , 



where the regions indicated in (5.1) are implied. The 
chemical potential serves as a Lagrange multiplier to keep 
the concentration conserved, so we can find /i by integrat- 
ing the concentration field over the size of the system and 
setting it equal to the average concentration (j)av- 



2L 



(t>siy)dy = (fia 
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We want the volume fraction X/3 of the background phase 
with concentration = —1 to be a;^ = {2L — X)/2L. Us- 
ing the lever rule 



and the equilibrium concentrations 
find that 



= 1, 



= — 1 we 



-1 



A 
T 



Doing the integral over the stationary concentration and 
keeping only the first order corrections in exp(— A) for 
A 1, we find that 



1 

M---e- 



(5.2) 



so that /i ^ as the system size _L — > cx). The depen- 
dence of /i on A will be important to our understanding 
of the physics in Section VI A below. 

Next we want to solve the full eigenvalue equation 



(I 



11) for the lamellar domain. Any perturbation of the 
domain can be written in terms of two linearly indepen- 
dent perturbation modes: either the two interfaces can 



fluctuate in phase with each other to form a "zig-zag" 
mode, or they can fluctuate out of phase in a "varicose" 
or "peristaltic" mode. These modes are pictured in Fig. 
^ and Fig. |3[ Since we are interested in calculating the 
eigenvalues near the marginally stable mode with uj — Q 
at fc = 0, we take the perturbed concentration field for 
the ziz-zag and varicose modes to be, respectively, 

0, = isech2(y + A/2) - isGch'(y - A/2) , (5.3) 



= -sech2(y + A/2) + -sech'(y - A/2) . 



The variational theorem (3.16) gives the eigenvalues for 
these two modes in the absence of any hydrodynamic 
effects. However, we are interested in the effect of the 
shear flow and of the fluid flow induced by gradients in 
the concentration. We cannot use the perturbation the- 
ory approach used in Section |^ because the varicose 
mode i s not a solution to any "unperturbed" operator 
in Eq. ( |3.1l|) (note that the zig-zag mode is the tran sla- 
tion mode, — 4>'g, and is an exact solution to ( 3.12| ) at 
k — 0). Instead, we adopt a "tight-binding" approxima- 
tion that will allow us to solve the full problem. 





FIG. 3. Varicose mode 



To implement this approach, we consider the two per- 
turbation modes above to be tw o basis states, and rewrite 
the eigenvalue equation ( 3.11 ) as a two- by- two matrix 
equation in this basis. We use as the right-hand ba- 
sis state and the conjugate function as the left-hand 



state. We insert our two trial functions (5.3) and (5.4) 
for (f) into the eigenvalue equation, multiply on the left by 
the corresponding 0, and integrate over all y. In vector 
notation, we have 



a(sech2(y + A/2) - sech^{y - A/2))/2 A 
6(sech^(y + A/2) -f sech^{y - A/2))/2 ) 

a(t>z 



(5.5) 



as our trial function, where a and h are the amplitudes 
of the two modes. Substituting into Eq. (3.11) gives the 
matrix equation 
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+ {(j>v,ik'yy(f>z) 



{4>z,vl(j)'^) + {4>z,ikjy(t)v) 



(5.6) 



Here we have used the definition rcf) ~ (p. The su- 
perscript on Vy indicates with which perturbation mode 
the velocity field corresponds, so that Vy is the velocity 
induced by the zig-zag mode and Vy the velocity indu ced 
by the varicose mode. On the left-hand side of ( |5.6| ) we 
use the orthogonality properties 

{4>z,4>v) = {4>vAz) = . 

These also apply to the diffusive terms on the right-hand 
side; this procedure thus ensures that in the absence of 
any flow effects we obtain the same eigen values lu as we 
would from the variational theorem ( |3.16 ). We can now 
solve (5.6) for the stability eigenvalues. Note that all cal- 
culations presented below are carried out to the lowest 
possible order in k. 



VI. LAMELLAR DOMAIN RESULTS 



A. Without Shear Flow 



We consider first the solution of (5.6) in the absence 
of the external shear flow (7 = 0). The only possible off- 
diagonal terms are the ones involving Vy. We begin by 
calculating the necessary integrals that form the matrix 
elements. 

Using Eq. ( 3.19| ) and expanding for small k as in Sec- 



tion 



mode is 



IV, we find the conjugate function for the zig-zag 



(j^ziv) = -lncosh(y+A/2)+lncosh(y-A/2)+fcA?/+0(fc^) 
The normalization integral is then 



(0,,0,)= / dy{(-lncosh(2/ + A/2) 

^ — 00 

+ lncosh(y - A/2) + k\y + 0(fc2)) 
X (sech2(y + A/2) + sech^{y - A/2)) / 2} 
= 2A - 2 - fcA^ + 0(fc2) . (6.1) 

Note that the second term on the right hand side of the 
above is negligible for sufficiently large A, but not when 
A is of the order of a few correlation lengths. Since it is 
reasonable to consider the case of A being a few times ^ 
(recall ^ = 1), we consider exp(— A) to be a small param- 
eter in the calculation, but not 1/A, so that we retain 
terms lik e the additive 2 in (3T). Next, substituting 0^ 
into Eq. (3.1C) for the velocity field, we find by expanding 
for small k as before, 



2 j + 0(fc4) . (6.2) 
s, there is only one term in F(f>z- 



Finally, since (pz 

F(f>z = ^k^{sech^{y + A/2) - sech^{y - A/2)) . (6.3) 

Now we turn to the varicose mode. The conjugate 
function for this mode, expanded for small k, is 



(pviy) = -J — lncosh(j/ + A/2) — lncosh(?; — A/2) 



+k [y' + ^ + Y^)+ Oik') . 

This leads to a normalization integral of 



(</>«, ^«) 



2A-2-h41n2-Hfc A' 



TT 
"3 



We note that the normalization goes to infinity as fc ^ 0. 
This is the mathematical manifestation of the fact that 
the varicose mode is not allowed at A: = because it does 
not conserve mass. For any nonzero k however there is 
no problem. The velocity field for the varicose mode is 
given by 

v^iy) = ^{2X- 3) e-^'y - + 0{k^e-^\ k^) . 

(6.4) 

In this expression we have not included terms of 
0(fc^e~^^). These terms will be negligible compared to 
terms of 0{k^) for k > e~^^ (at such small k, of course, 
the linear term in k will dominate over any k^ terms). 
We will see below that this condition is met for the k 
values of most interest. Finally, for the varicose mode 

Fcl)^ = ifc2(sech2(y + A/2) + sech2(y ~ A/2)) 

-3sech^(t/ + A/2)sech2(y - A/2) , (6.5) 

so that Fcjjy includes an overlap term between the two 
interfaces. 

It is fairly simple to show by strai ghtf orward inte- 
gration that the off-diagonal terms in (5.£) vanish (for 
7 = 0): 



^^'^y't's) ~ vvuj "y 



'.,<0l.) = 0. 



This reduces the matrix equation to 
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,(t>z)uJ \ ( ^ 

{(^,,(j>,)Lu J \b 







(6.6) 



so that we can solve for each eigenvalue separately 



(6.7) 



for each mode. Using the expressions given above we per- 
form the remaining integrals to obtain w for each mode. 
For the zig-zag mode we find 



37? 6(A- 1) 



1 



fcA^ 
2A - 2 



where /(A) is the function 

^ ^ ^2(A)-AJi(A)-^o 
^^^>- A^l ■ 

Here 5q is the definite integral 



6o ^ dyy^sech^ylncoshy = 1.70681, 
and the functions 6i and 62 are the overlap integrals 

/CO 
dyysech^y\ncosh{y — A) , 
-00 

/oo 
c?2/y^sech^?;lncosh(?/ — A) . 
-00 

These are integrated numerically using Mathematica. We 
see from Eq. (6.8) that the zig-zag mode is stable for small 



k. The terms involving the dimensionless viscosity 77 are 
due to the flow field induced by the perturbations in th e 
concentration cj), and come from the Vy term in (3.11). 
Depending on the value of ^, either the hydrodynamic 
terms or the diffusive terms will dominate. We find that 
the stability eigenvalue for the varicose mode is 



„, 2A 12fcAe-2^ /2 



V 
or] 



fc3 

12 



(6.9) 



The varicose mode is thus unstable for sufficiently small 
wavelengths! The eigenvalues for the two modes are plot- 
ted as functions of k in Fig. |^. Here we take A = 6, so 
that e~^^ is small (as we assumed above) and k > e~^^ 
for most of the range in the graph, as discussed above. 
We take the dimensionless viscosity to be = .1, which 



is a typical value for critical binary fluids. However the 
overall shape of the dispersion relations remains similar 
for other values of A and fj. 

The instability of the varicose mode may seem unintu- 
itive. We first note that it is unstable only for sufficiently 
small fc, and is stabilized at larger wavenumbers by the k^ 
curv ature term, the same term that was obtained in Eq. 
(4.2) for a system with a single interface. Second, the in- 
stability is exponentially small in the separation between 
the interfaces, A. This is thus a very weak instability. It 
is due to a coarsening effect (essentially Ostwald ripen- 
ing) in which thin regions of the middle phase shrink in 
favor of fatter regions. Recall that the chemical potential 
H ~ e""**. If A decreases in a region, /i increases, so the 
chemical potential is higher in the neck regions than in 
the bulges. This drives a flux from the necks towards the 
bulges (see Fig. We can understand the lowest order 
diffusive effect as follows. First note that the lowest order 



diffusion term in Eq. (x£), with units, is 



oj l6D-e-^^^^ 

As before, we can express the velocity of the interface uiA 
as 

wA(j)e ~ DV(j) - Dk6(j) 

since the concentration gradient is along the x-direction. 
The excess concentration added (subtracted) in the bulk 
regions of the necks (bulges) is essentially 



sech^y/^ 



-2A/C 



so that 



k 



This implies that a large sheet of one phase immersed 
in the other will break up into cylinders via this insta- 
bility. Note that this is not the Rayleigh instability of a 
long fluid cylinder, in which the cylinder is unstable to- 
wards long wavelength, axisymmetric fluctuations. That 
is a hydrodynamic instability that occurs for a three- 
dimensional cylindrical interface because the curvature 
at the necks is higher than at the bulges. In this two- 
dimensional perturbation mode, the curvature at the 
necks and bulges is of the same magnitude (the extra 
dimension out of the plane of say Fig. |5| does not exist), 
and so there is no curvature-driven instability. The cur- 
vature effect is stabilizing, and it is the thermodynamic 
force driving phase separation that causes the instability. 
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FIG. 4. Dispersion relation for 7 = 0, with A = 6 and rj = .1 



X 




FIG. 5. Diffusional instability of the varicose mode 



B. With Shear Flow 



Next we consider what happens when we include the These two off-diagonal elements are found to be 
external shear flow. Physically, the shear flow tends to 
mix the two modes since the top interface travels in an 
opposite direction to the bottom interface. We might 
then expect that at some shear rate, the two perturba- 
tion modes lose their distinguishing features. 

To calculate the eigenvalues we only need to calculate 
the matrix elements involving the shear. It is straight- 
forward to show that the operator ikjy is off-diagonal in 
the basis of our two perturbation modes, i.e. 



-0(fc3), 



j,ik'yy(l)^^) = -2i^\ + ik-^{\ - 2Aln2 - 5i + A^) 



^fc^7(y + ^)+0(fc3). (6. 



The stability eigenvalues are now found by diagonalizing Eq. ( |5.6| ), which means solving the secular equation 

= 0. (6. 



10) 



11) 



w^-w^ {4)z,ik'^y(t>v)/{4>z,(t)z) 

{4>y,ikjy(l)z)/{(i>v,4'v) uj.y-uj 



12) 



Solving for u gives 



1 1 / 2 

w±(fc) = -{oJz{k) + ujy{k)) ± -y {<^z{k) - uJv{k)) - 72s(fc) , 



where s{k) is given by 

k^X{X^-X-Si) 



5(fc) = 



fc3 



A- 1 



2(A-1) 



'^^ + A^) (A^-A-<5i)-A^-yA2 



+ 0(fc4). 



(6.13) 



(6.14) 
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Some examples of the two curves Re{uj±{k)) are shown 
in Fig. ^ The spacing between the two interfaces is 
A = 6£ an d we have taken fj — .1. For 7 = it is clear 
that (|03|) reduces to our previous results, with uj+ = lo^ 
and uj^ ~ Uv Fig. ^ shows lj_ for three different shear 
rates (it turns out that the curves for uj^ for these same 
shear rates are nearly indistinguishable, so they are plot- 



ted as one curve in Fig. |6|). We see that at low shear 
rates the unstable mode still exists, but the window of 
wavenumbers over which Re{LoJ) < becomes smaller 
as 7 increases. Above some critical shear rate 7c, the 
previously unstable mode becomes stable for all k. The 
shear flow thus completely stabilizes the varicose mode, 
by mixing it with the stable zig-zag mode. 



0.1 




FIG. 6. Dispersion relations: 
(dash-dot) 



0.02 



UJ+ for 7 = .04 (solid line); 100a;_ for 7 = .04 (dashed), 7 = > = .07225 (dotted), and 7 : 



We can easily solve for t he cri tical shear rate First 
note that the first term in (6.13) is positive, because the 
negative terms in ujy are exponentially small in A. As 7 is 
increased, the square root term in (S.13) becomes smaller. 
The effect is that the value of k below which < be- 
comes smaller with increasing shear; the domain is only 
unstable to longer and longer wavelength perturbations 
as the shear rate is increased. For a given shear rate, 
u!^ > for all k > kc where kc satisfies a>_(fcc) = 0: 



7c 



4(A~ l)c-^^(8?? + 4A(2A-3)) 
3^2A(A2 - A + 5i(A)) 



(6.17) 



UJz{kc) +cj„(fcc) = \/ {LOz{kc) - oJvikc))^ - 7^s(fcc) ■ 

(6.15) 

The unstable mode becomes stable for all wavenumbers 
k when kr — + 0. To find the critical shear rate, we first 



solve ( 6.15 ) for ^{kc): 

.2/, X _ -4:UJz{kc)uJv{kc) 



s{kc) 



(6.16) 



Taking the limit kc ^ in Eq. ( |6.16 ) gives the critical 
shear rate for complete stabilization: 



For the specific values A = 6 and 77 — .1, one finds 
7c — .07225 as indicated in Fig. ^. 

The critical shear rate is graphed as a function of 77 
and A in Fig. |^ and Fig. ^. We note the lamella is sta- 
ble for all 7 > 7c. We see that 7c is an algebraically 
decreasing function of fj and an exponentially decr easing 
function of A. Recall that fj — ADri/3a£,, so that ( 6.17 ) 
tells us that as the viscosity increases, or alternatively as 
the surface tension decreases, the easier it is for the shear 
flow to mix the two modes before the unstable perturba- 
tion has a chance to grow. We can also invert ( 6.17 ) to 
obtain the critical width Ac above which the lamella is 
stable for a given shear rate 7. As we see from Fig. ||, 
given a shear rate 7, at values of A lying below the curve 
the lamella is unstable to the varicose coarsening mode 
whereas for values of A above the curve the lamella is 
stable and will no longer coarsen. This simple system of 
a single lamellar domain thus exhibits the well-known ex- 
perimental observation that the shear flow tends to halt 
the phase separation process. 
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FIG. 7. Critical shear rate 7c (fy) for A = 6 
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FIG. 8. Critical shear rate 7c(A) for r] — .1 



VII. DISCUSSION 



We have seen that in the case of an isolated lamellar 
domain, shear flow has the effect of mixing the zig-zag 
and varicose modes so that they both become stable. Es- 
sentially, the flow eliminates the special phase relation- 
ship between the two interfaces necessary for the vari- 
cose mode to exist. The physics of this mode is that 
thin regions evaporate in favor of thick regions, but in 
the presence of shear thin and thick regions do not ex- 
ist long enough for this diffusion to take place since the 
fluctuations are being carried downstream. 

We would expect that a similar mechanism would ap- 
ply to a large stack (along the y-direction) of lamellar 
domains. Although the stability eigenvalues have not 
been calculated for this case, the effects seen in the sin- 



gle lamellar domain should apply. Coarsening in the y 
direction in a stack of lamellae is also dependent on thin- 
ner regions evaporating, their atoms diffusing across the 
intervening phase to a thicker region. From we ex- 
pect this coarsening instability to also have a rate that is 
exponentially small in A. When one considers sinusoidal 
perturbations of the layers in a shear flow, once again 
the phase relations between interfaces will be constantly 
changing. As A increases, the atoms must diffuse far- 
ther across a layer for the pattern to coarsen, but they 
must be able to do so before they are swept downstream 
by the shear flow to a new x position where the diffu- 
sion is no longer favored. We might anticipate then that 
in a general two-dimensional system with many lamellar 
domains, for any given shear rate 7 there is an upper 
limit Ac to the layer spacing for which the coarsening 
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instability is still present. The shear flow destroys the 
correlations between interfaces necessary for the coars- 
ening instability to operate, leading to a dynamic steady 
state. The strength of the shear flow would determine 
the typical lamellar width Ac (7) present in the system at 
steady state. 

This behavior is qualitatively similar to that seen in the 
fully three dimensional "string" phase in shear flow. We 
do not expect quantitative agreement, however, because 
the stability analysis of the lamellar domain considered 
here is strongly dependent on the dimensionality. The 
instability of a long cylinder is much stronger than the 
weak exponential (2D) instability found here. For the 
case of a viscous cylinder of fluid immersed in another 
viscous liquid, the hydrodynamic instability correspond- 
ing to a varicose perturbation has a dispersion relation 
that behaves as [|l6| 

uj :^f{ka) 

Zrja 

where a is the radius of the cylinder. Thus, we might 
expect more dramatic effects in this case. 

In summary, we have shown that a long extended do- 
main in the two-phase state of a two-dimensional, phase- 
separating binary fluid can be stabilized by an applied 
shear flow. There is a critical shear rate below which 
the extended domain is unstable towards long wavelength 
fluctuations and above which we predict complete stabi- 
lization. This is in qualitative agreement with experi- 
ments on dynamic steady states in phase-separating flu- 
ids under shear flow, however the mechanisms operative 
here are different due to the reduced dimensionality. We 
intend to report results of a similar calculation for a long 
cylindrical domain under flow in the future. 
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